Non-equilibrium electronic structure of interacting single-molecule nanojunctions: 
vertex corrections and polarization effects for the electron-vibron coupling. 
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We consider the interaction between electrons and molecular vibrations in the context of electronic 
transport in nanoscale devices. We present a method based on non-equilibrium Green's functions to 
calculate both equilibrium and non-equilibrium electronic properties of a single-molecule junction in 
the presence of electron-vibron interactions. We apply our method to a model system consisting of 
a single electronic level coupled to a single vibration mode in the molecule, which is in contact with 
two electron reservoirs. Higher-order diagrams beyond the usual self-consistent Born approximation 
(SCBA) are included in the calculations. In this paper we consider the effects of the double-exchange 
diagram and the diagram in which the vibron propagator is renormalized by one electron-hole bubble. 
We study in detail the effects of the first- and second-order diagrams on the spectral functions for 
a large set of parameters and for different transport regimes (resonant and off-resonant cases) , both 
at equilibrium and in the presence of a finite applied bias. We also study the linear response (linear 
conductance) of the nanoj unction for all the different regimes. We find that it is indeed necessary 
to go beyond the SCBA in order to obtain correct results for a wide range of parameters. 

PACS numbers: PACS numbers: 71.38.-k, 73.40.Gk, 85.65.-(-h, 73.63.-b 
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I. INTRODUCTION 

Single-molecule electronics has shown significant 
progress in the recent years. A variety of interesting 
effects have been observed in the transport properties 
of single (or a few) conjugated organic molecules in- 
cluding rectification, negative differential resistance, and 
switchingi"— . In these quasi-one dimensional systems, 
which present well delocalised 7r-electrons, the electronic 
current fiowing through the quite fiexible backbone of the 
molecule affects the ground state properties of both elec- 
tronic and mechanical degrees of freedom of the molecule. 

The importance of inelastic effects in the transport 
properties has been demonstrated in several ground- 
breaking experiments^"^ ^; these effects are related to the 
interaction between electron and mechanical degrees of 
freedom of the molecule. 

The interaction between an injected charge carrier 
(electron or hole) and the mechanical degrees of freedom 
(phonon, vibron) in the molecular junctions is impor- 
tant in order to understand energy transfer, heating and 
dissipation in the nanojunction^^. The electron-vibron 
interaction is also at the heart of inelastic electron tun- 
neling spectroscopy (lETS). lETS is a solid-state-based 
spectroscopy which gives information about the vibration 
modes (vibrons) of the molecules in the nanoj unctioni^. 
It is now possible to measure such vibrational spectra for 
single molecules by using scanning tunneling microscopy 
(STM)aii to build lETS map&i^., or by using other elec- 
tromigrated junctions or mechanically controlled break 
junctionaii^iiaiii. 

At low applied bias (typically 100-400 meV) the lETS 
show features (peaks, dips, or peak-dip-like lineshape) 
which have been attributed to selective excitation of spe- 
cific vibration modes of the molecule. The position in en- 



ergy (bias) of the features correspond approximately to 
the frequency of the vibration, as given by other spectro- 
scopic data (IR, Raman) obtained on the same molecules 
in a different environment. 

There have been many theoretical investigations focus- 
ing on the effects of electron-vibron coupling in molecu- 
lar and atomic scale wire^i^"— . Most of them have fo- 
cused on the interpretation of the features in lETS. How- 
ever, most of these studies have been performed by us- 
ing the lowest-order expansion possible for treating the 
effects of the electron-vibron interaction (i.e using the 
so-called self-consistent Born approximation SCBA). In 
the language of many-body perturbation theory, it cor- 
responds to a self-consistent Hartree-Fock calculation for 
the electron-vibron coupling. 

However, in analogy to what is obtained at the 
Hartree-Fock level for interacting electrons, there are 
good reasons to believe that this approximation is not 
enough to correctly describe the physics of the electron- 
vibron interacting system, especially beyond the weak 
electron-vibron coupling regime. For example, the limits 
of SCBA have already been investigated in Ref.jsO] but 
without introducing remedies to go beyond SCBA. 

In this paper, we examine this using a 
true non-equilibrium Green's-function (NEGF) 
technique20.>25^>28^-^^d4-i6 ^hich allows us to 
study all the different transport regimes in the presence 
of electron-vibron interaction. Following the spirit 
of many-body perturbation theory and Feynman dia- 
grammatics, we go beyond the commonly-used SCBA 
approximation by introducing higher-order diagrams for 
the electron-vibron interaction. 

We study the simplest possible model system which 
nonetheless contains the relevant physics of the trans- 
port properties of the molecular junction2S,i22. Further- 
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more, because of the uncertainty of the exact geometry 
of the single-molecule junction in the experiments, there 
is a corresponding uncertainty about how to model the 
coupling between the molecule and the electrodes and 
correspondly for the potential drops at each molecule- 
electrode contacts. Hence we take the quantitites char- 
acterizing the potential drops at the contacts as phe- 
nomenological parameters^^i^. 



We concentrate in this paper on the electronic prop- 
erties of the molecular junction in both equilibrium 
and non-equilibrium conditions as well as on the linear- 
response properties of the junction (prior to considering 
the full non-equilibrium transport properties in a forth- 
coming paper) . Such properties are given by the density 
of electronic states and represented by the spectral func- 
tions, which are at the very heart of all physical proper- 
ties of the system, such as the charge density, the current 
density, the total energy, etc. 




FIG. 1: Schematic representation of the single site single 
mode (SSSM) model. The single electronic level eo is broad- 
ened by the coupling to the leads; the electronic transport is 
resonant if the Fermi levels fiL,R are very close to this broad- 
ened level, and off-resonant (as shown) otherwise. The system 
is shown in equilibrium, with hl = ftR, non-equilibrium ef- 
fects are studied by moving the Fermi levels. 



Spectral functions are most closely related to photoe- 
mission and adsorption spectroscopies. To our knowl- 
edge such experiments have not yet been performed on 
single-molecule junctions, though photoemission spectra 
have been measured on quasi one-dimensional supported 
atomic scale metallic wires (see for example Ref.jg^l 
showing interesting results on one-dimensional collective 
electronic excitations) . 



The paper is structured as follows. We start with 
a description of our model system in Section III Al and 
a discussion of the relevant underlying theory of non- 
equilibrium Green's functions in section IIIBl Our cal- 
culated spectral functions are presented in Section lllli 
where we consider first the equilibrium case (Section 
IIII Ap and then the non-equilibrium case (Section IIII Bp 
at the Hartree-Fock level. We discuss especially the ef- 
fects of including or not the Hartree diagram in the cal- 
culation. We also compare NEGF calculations with re- 
sults obtained from inelastic scattering techniquesi^i^iH 
for equivalent model systems in Section fill CI We show 
that it is indeed necessary to go beyond SCBA to obtain 
correct results for the relevant range of electron-vibron 
coupling. We then present the effects of the second- 
order diagrams in the spectral functions in Section ITlI D I 
The second-order diagrams correspond to two classes of 
process; the first is related to vertex corrections of the 
SCBA calculation and the second to polarisation effects 
(i.e. partial dressing of the vibron propagator by the 
electron- hole bubble diagram) . In Section IIII El we dis- 
cuss the effects of different levels of approximation for the 
electron-vibron coupling on the linear conductance of the 
molecular junctions. Throughout the paper, we will use 
the term vibron to define a quantum of vibration of a 
mechanical degree of freedom. 



II. MODEL 
A. Hamiltonian 

Our model is based on a system with an interact- 
ing central region connected to two non-interacting elec- 
trodes (see figure [Ij. 

The total Hamiltonian for the system is given by 

^^total = -H^e + ffvib + ^^e-vib, (1) 

where He,i7vib and i?e-vib represent the electronic, vi- 
bron, and electron-vibron coupling parts of the Hamil- 
tonian respectively. The electronic part of the Hamil- 
tonian is broken into sections describing the left (right) 
electrode {Hn), the central interacting part Hq and 
the potentials coupling the central region to the left and 
right electrodes Vlc + VcR respectively: 

H, = Hl+H^ + Hr + Vlc + Vcr. (2) 

The Hamiltonians for the electrodes are given by 

Hl + Hr= ^ eaclca, (3) 

a=L.,R 

where cj^ (ca) creates (annihilates) a non-interacting 
electron with energy e on electrode a. The electronic 
Hamiltonian for the central region and the coupling po- 
tentials are given by 

^c = ^c({4},K}), (4) 

VlC + VcR^ J2 Vo.ncidn + V*^dica, (5) 
n,a.—L,R 
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where the interacting electrons in the central region are 
created (annihilated) on electronic level n by the opera- 
tors djj (c?n). 

In principle, our Hamiltonian for the central region 
may contain electron-electron interactions, and is built 
from a complete set of single-electron creation and anni- 
hilation operators. In our current work we do not include 
any electron-electron interactions, and so the electronic 
part of the total Hamiltonian for the central region He 
becomes 

Hc^^endidn- (6) 

n 

Meanwhile, the vibron part of the Hamiltonian is rep- 
resented by 

^^vib = ^ huxalax, (7) 

A 

where a\ creates (annihilates) a vibron in vibron mode 
A with frequency uj\. The electron- vibron coupling term 
is taken to be linear in the vibron displacement and its 
most general expression is then given byi^ 

-ffe-vib = ^ IXnrnial + ax)dl^dm, (8) 
A,n,m 

where "fxnm is the coupling constant for exciting the vi- 
bron mode A by electronic transition between the elec- 
tronic levels n and m. 

We concentrate on the simplest version of the Hamil- 
tonian of the central part: the single-site single-mode 
SSSM model, in which one considers just one electron 
level coupled to one vibration mode. Despite the sim- 
plicity of this model, remarkably it not only contains all 
the physics we require but also ensures that we isolate the 
properties we are interested in without the complication 
of the added electronic levels. The reasons for this are 
as follows. Firstly, when the Fermi energy of the leads is 
pinned around the midgap, then at low and intermedi- 
ate biases one of the frontier orbitals (either HOMO or 
LUMO) dominates the transport properties. Secondly 
in conjugated organic molecules (mostly used in single- 
molecule junction experiments), it is known that the op- 
tically relevant vibration modes are mostly coupled to 
either the HOMO or LUMO levels^. We can therefore 
be confident that our model contains the relevant physics. 

The total Hamiltonian for the central region thus be- 
comes 

He = £od^d + huJoa'''a + "fo{a''' + a)d^d, (9) 

where we now have just one electronic level e„ — >■ Sq and 
just one vibron mode ax a, coupled via the electron- 
vibron coupling constant 70 . The expression for the lead- 
central- region coupling (equation ([5])) also simplifies to 
become 

Vlc + Vcr^ tocicU + d^Ca): (10) 

a=L,R 

where we have replaced the coupling potentials Van with 
hopping integrals ioa- 



B. Non-equilibrium electron Green's functions and 
electron-vibron self-energies 

Non-equilibrium Green's functions (NEGF) within the 
Keldysh formalism^""— represent an extremely useful 
tool for studying the non-equilibrium properties of many- 
particle systems. The Green's functions are functions of 
two space-time coordinates, and are obviously more com- 
plicated than the one-particle density which is the main 
ingredient of density-functional-based theories. One of 
the great advantage of NEGF techniques is that one can 
improve the calculations in a systematic way by taking 
into account specific physical processes (represented by 
Feynman diagrams) which is what we do in this paper for 
the electron-vibron interaction. The Green's functions 
provide us directly with all expectation values of one- 
body operators (such as the density and the current), 
and also the total energy, the response functions, spec- 
tral functions, etc. 

In Appendix 13 we provide more details about NEGF 
and how to obtained the electron-vibron self-energies 
from a Feynman diagrammatic expansion of the electron- 
vibron interaction. We now briefly describe how we apply 
the NEGF formalism to the SSSM model. 

The Green's functions are calculated via Dyson-like 
equations. For the retarded and advanced Green's func- 
tions C'" these are 

G"" - Qc" + 9c"^'''''G'-^\ (11) 

where g^"' is the non-interacting Green's function for the 
isolated central region. 

For the greater and lesser Green's functions, 
we use a quantum kinetic equation of the form 

G>^< = (1 + G''I]'')5>'<(1 + S°G") + G''S>'<G°. (12) 

Here S^, (a; = r, a, >, <) is a total self-energy consisting 
of a sum of the self-energies from the constituent parts 
of the system: 

I]- = Ei + S5,-|-Sf„,. (13) 

R self-energies arising from the non- 

interacting leads a = L,R and as such are simple to 
calculate: 

S^^t^gSJc.), (14) 
= (KT, (15) 

E>-2z3mK(c.)] (l-/„(c.)), (16) 
I]<--2*5mK(^)] /„(..), (17) 

where is the Fermi-Dirac distribution for lead a, with 
Fermi level = + VaS^V and temperature Tq. The 
fraction of potential drop at the left contact is 77^ = ±77^ 
and i]R — =f(1 — rjv) at the right contact^^, hence rji — 
VR = 6^ is indeed the applied bias, and rjv S [0, 1]. 

The component of the retarded Green's function for 
the isolated (non-interacting) lead a corresponding to the 
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FIG. 2: The (a) Fock and (b) Hartree diagrams 



(a) SjJ^^iij 




(b) 




FIG. 3: The (a) double exchange DX and (b) vibron prop- 
agator dressed by the e-h bubble diagrams (dressed phonon, 
DPH) 



sites (or energy levels, depending on the representation 
used to the electrodes) connected to the central region is 
given by gl^. 

In this paper, we have chosen a simple model, a semi- 
infinite tight-binding chain with on-site energy and 
nearest-neighbour hopping integral /3q. . This model gives 
a semi-elliptic density of states of the terminal lead 
site connected to the central region, and so each lead's 
Green's function becomes 



9l(Sy>) = exp(ifca(a;))//3o 



(18) 



with a; = £q + 2/3^ cos ka{uj). We have chosen this model 
because it is one of the most simple, although in prin- 
ciple and in practice there are no limitations for taking 
any other more complicated or more realistic models for 
the lead, such as Bethe lattices with z-coordination, or a 
nanotip supported by a semi-infinite surface as shown in 
Ref |6J) since all their electronic properties are wrapped 
up in the lead self-energies E| r{^)- 

The self-energy for the interacting central region, Eint 
is somewhat more complicated. It consists of the sum 
of the self-energies due to interactions between the elec- 
trons and to interactions between the electrons and the 
quantum vibration modes (vibrons). In this paper we 
consider only the coupling between each electron and a 
single vibron of the central region (the molecule) , hence 
the electron- vibron self-energy Eg-vib- 

For the current work, it is necessary to calculate sev- 
eral types of self energy. Firstly we have the Fock- 
like self-energy 'E'^^^/^^^{uj), which is a function of en- 
ergy, and is represented in diagrammatic form by Fig- 
ure mja). We also have the Hartree- like self-energy 



which is independent of energy, given by Figure 
[2Ub). Calculations using only the Hartree and Fock di- 
agrams, and performed in a self-consistent way are usu- 
ally referred to as the self-consistent Born approximation 

( S C B A ) 20.24. 26. 27. 30-34. 38.44-46 

However, as explained in the introduction, we also 
want to go beyond the SCBA, and thus we will also cal- 
culate two further self-energies that include two-vibron 



processes. The first of these is the double-exchange self- 
energy S^^jjj given by Figure ^a) . In the many-body 
language, it is part of the vertex correction to the Fock 
diagram. The second is given by Figure and cor- 
responds to the dressed vibron (or GW-like) self-energy 
^e^vib' "which consists of the vibron propagator renormal- 
ized by a single electron- hole bubble (the polarization). 
This is why we refer to the effects of ^^XR polarization 
effects in the following. The details of how we implement 
these self-energies are given in appendix 1X1 



C. Physical properties 

Once we have calculated all the different Green's func- 
tions, any of the physical properties of the system, such 
as the electron density, the electronic current denstiy, the 
total energy, the current noise, the heat transfer, etc. can 
be calculated. 

For example, the electronic current passing through 
the a contact is given by 



{ci{t)d^{t))-v:jdi{t)ca{t)), (19) 



i.e the first term describes the transfer of an electron 
from the interacting region to electrode a, while the sec- 
ond transfers an electron from electrode a to the central 
region. 

We can then express this in terms of Green's functions 
and derive an expression for the expectation value of the 
current2fi: 

/a - 11 ^ Tr{S<(c.) G>(c.) - S>(c.) G<(c.)} 

= ^Tr{/„(a.)r„(..)[G'-(^)-G''(.;)] 
+ r„(w) G<(w)} . 

(20) 

All physical properties may be expressed in terms of 
the spectral function A{u!) which is at the heart of this 
paper. The spectral function is related to the imaginary 
part of the retarded or advanced electron Green's func- 
tions, as 



A{lu) = -3m[G'"(a;)] = +9to[G''(w)] 



(21) 



For non-interacting systems, it is simply proportional to 
the density of electronic states n{uj) — A{u!) / tt. For inter- 
acting systems, it gives information about the excitations 
(electron or hole) of the system. 

Furthermore, when the system is at equilibrium (/^ = 
fn = /°'i), there are some relationships between the 
lesser (greater) and the advanced and retarded Green's 
functions: 



G<^'='i(w) ^ -r'^ioj) (G'''°'i(w) - C'^'^iuj)) 



(22) 
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and 

G>'^^{u:) = -iniuj) 1) (G^'^'^ioj) C'-^^iuj)) 

= 2i(r»-l)AH. ^ > 

These relationships are at the centre of the fluctuation- 
dissipation theorem for equilibrium, which also be recast 
as a relationship betweeen the greater and lesser Green's 
functions: 

G>^'='5(cj) = ^(,i'^-^^o)/kT G«,oq(^) (24) 

for statistical averages at finite temperature in the grand 
canonical ensemble. This equation is related to the 
Kubo-Martin-Schwinger boundary conditions . 

For non-equilibrium conditions, there is no unique 
Fermi level at finite bias (or no unique temperature if 
Ti 7^ Tr) in the whole system, and the relationships 
given by Eqs. (|22j|24p no longer hold. This is an impor- 
tant feature of the non-equilibrium formalism for which 
conventional equilibrium statistics need to be reformu- 
lated. 



D. Computational aspects 

The calculations start by constructing the non- 
interacting Green's functions of the entire system Gg 
with X being any three of the possible Green's functions 
X = {r,a, <,>,t,t}. The other Green's functions are 
obtained by using the relationships between them as 
shown in Appendix lA II For example, the retarded 
Green's function is given by Gq{uj) ~ [ g5(w)^^ — I]£(a;) — 

One then calculates three "Keldysh components" for 
the self-energies corresponding to any of the diagrams 
★ = {F,DX or DPH}: E*;< , S*;> and^S*;*,. The 
Hartree diagram has only one component as shown 

above. 

The advanced and retarded self-energies S^ ,^.^' ^ are 
then obtained by simple algebra using the relationships 
between the different self-energies as explained in Ap- 
pendix lA 11 

The new Green's functions are then calculated by us- 
ing the Dyson equations for G^''^ 

G^''\u) = [u^-eo-K':a^)r\ (25) 

and the quantum kinetic equations for G^'^ 

G^{u;)^G^{uj)j:ia^)G'^i^)^ (26) 
where the total self-energies are given by 

any* 

The new self-energies are then recalculated and the pro- 
cess re-iterates until full self-consistency is achieved. 



Actually at each iteration of the calculations, we use a 
simple mixing scheme of the self-energies obtained at the 
present iteration and at the previous iteration. This mix- 
ing scheme permits us to achieve full self-consistency in a 
maybe slightly longer but more stable iterative process. 

Note finally that by using Eq. [211 instead of the more 
general formulation given by Eq. [T21 we assume that, 
after switching on the interactions, there are no bound 
states in the system (i.e. there are no interaction-induced 
electron states located outside the spectral supports of 
the left and right leads)^^'^"^ which is indeed the case. 



III. RESULTS: SPECTRAL FUNCTIONS 

In this section, we present results for the spectral func- 
tions A{uj) for the different transport regimes and for 
different applied biases in the low vibron-temperature 
regime. 

We divide the calculations into four types. Firstly, we 
calculate the spectral functions at equilibrium for two 
transport regimes. The first of these is where Eq 3> fJ-^^ 
or So <C known as the off- resonant regime as in order 
to create a current between the left and right leads one 
puts an electron in the empty (for electron transport, 
£o ^ m"*^) or full (for hole transport, eo ^ electronic 
level Eq. 

The second transport regime is when Eq ± linewidth ^ 
fi°'^, known as the resonant transport regime. The 
linewidth is the width of the peak in the spectral function 
which arises from the electronic coupling of the central 
region to the left and right leads. In this case the elec- 
tronic level £o is, at equilibrium, half-filled by an electron 
(and thus also half-filled with a hole). 

For each of these transport regimes, we calculate the 
spectral function both at equilibrium (applied bias V — 
0) and non-equilibrium at finite bias (^ > 0). 

For each of these four groups (resonant/off-resonant 
transport regimes, at/out of equilibrium) , a large amount 
of different NEGF calculations have been performed for 
different values of the relevant parameters and within 
different levels of approximation (Hartree Fock, Hartree 
Fock-|-second-order, partially or fully self-consistent cal- 
culations). In the rest of this section, we present only a 
limited and selected number of results which we found 
the most relevant for each case, and we analyse and com- 
pare in detail the effects of the different diagrams on the 
spectral functions of the system at and out of equilib- 
rium. 

We also, in section llll C[ compare perturbation expan- 
sion based calculations (NEGF) to a reference calcula- 
tion which is exact in term of electron- vibron coupling 
but which however is only valid for a specific transport 
regime. 
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FIG. 4: Equilibrium (zero bias) spectral functions A{ijj) for 
the off-resonant (a) electron and (b) hole transport regime. 
Calculations were performed with the Fock-like electron- 
vibron diagram (solid line) and with both the Hartree and 
Fock-like diagrams (symbols and dotted line). For the elec- 
tron transport regime the inclusion of the Hartree self-energy 
has no effect, but for the hole transport regime it shifts the 
entire spectral function to lower energies. The parameters are 
£o ~ -1-0. 5(— 0.5) for electron (hole) transport 70 — 0.21, a;o ~ 
0.3, toL,R = 0.15, 7? = 0.005 and hl = fiR ^ fJ-""^ = 0. 



tron excitations. The main peak corresponds to adding 
an electron in the single available level. This peak is 
broadened by the coupling to the leads, and its position 
in energy Eq is renormalised by the electron- vibron in- 
teraction, i.e. Eq is close to the static the polaron shift 

Eq^Eq- 'Jq/oJq. 

The lesser peaks in the electron-transport spectral 
function are vibron side-band peaks arising from reso- 
nance with n — 1,2,3,... excitations in the vibration 
mode. These peaks correspond to vibration excitation 
(vibron emission) only. At zero vibron temperature, 
these are the only available mechanisms for vibrational 
excitations. We note, and discuss further in section llil CI 
that these side-band peaks should occur at integer mul- 
tiples of ujQ away from the main peak, but that for both 
our Fock-only and Hartree-Fock SCBA calculations the 
peak-peak separation is slightly wider than this. 

In this regime, the Hartree self-energy is negligible, be- 
cause most of the spectral weight is above the Fermi level 
and S^^j^ oc dujA{uj) ^ 0. This implies that the po- 
laron shift is mainly due to the Fock-like self-energy. 

For the off-resonant hole-transport regime ( Fig- 
ureHl^b)), all the features in the spectral function occur at 
w < and therefore correspond to hole excitations. The 
vibron side band peaks are at lower frequencies than the 
main peak because they correspond to the emission of 
vibrons by holes rather than electrons. 

When we include just the Fock-like self-energy, the 
hole spectral function is symmetric (with respect to the 
equilibrium Fermi level fj,°'^) with the electron spectral 
function, as can be seen clearly in figure [5j Adding the 
Hartree self-energy, however, breaks this electron-hole 
symmetry. As most of the spectral weight is below /i'^'^, 
the expression for S^^jj^ given in equation Ea. (jA15|) re- 
duces to a constant (i.e. twice the polaron shift) 
as f |^iG^(a;) 1. As a result of this the whole spectral 
function is shifted to the left by this amount. 



A. At equilibrium 

We first consider the spectral functions at equilibrium, 
with no applied bias. Figure S] shows A{uj) for the off- 
resonant regime, and Figure [5] shows the resonant trans- 
port regime. 

In the equilibrium many-body language, the features 
in the spectral functions obtained at positive energies 
(w > 0, above the Fermi level) correspond to electron 
excitations, while the features at negative energies (a; < 
0, below the Fermi level) correspond to hole excitations. 



1. Off-resonant transport regime 

For the off-rcsonant electron transport regime (Fig- 
ure m^a)) all the features in the spectral function are 
above the Fermi level fj°'^ and hence correspond to elec- 



2. Resonant transport regime 

We now turn to the resonant transport regime, with 
the spectral function shown in Figure [5l Here our elec- 
tronic level eo is broadened by the coupling to the leads 
and is partially filled with electrons. We can see that cal- 
culations performed with only the Fock-like self-energy 
preserve the electron-hole symmetry. The spectral func- 
tion presents peaks located both at positive and negative 
energies, which correspond to the emission of vibrons by 
electrons or holes respectively. As in the previous sec- 
tion, the inclusion of the Hartree self-energy (equation 
(|A15I) ) breaks down the electron- hole symmetry and the 
features are shifted to lower energies with a correspond- 
ing modification of the spectral weights for each peak. 

We note here that the inclusion of the Hartree self- 
energy thus modifies so drastically the spectral function 
that it will also strongly affect the I-V characteristics 
of the junction in comparison to calculations performed 
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FIG. 5; Equilibrium spectral function A{lj) for the resonant 
transport regime. Calculations were performed with the Fock- 
like (solid line) and with both the Hartree and Fock-like di- 
agrams (dashed line). The parameters are identical to those 
used in Figure U except for the value of the electron level 
which is eo = for the resonant regime. 

with only the Fock self-energj«2S. 

B. Non-equilibrium spectral functions 

1. N on- equilibrium off-resonant transport 

In this section, we present results for the spectral func- 
tions A{io) for the off-resonant transport regimes, for dif- 
ferent applied biases, both with and without the Hartree 
contribution. We first present calculations for an asym- 
metric potential drop with ?7y = 1, i.e. fJ-L = V and 

m = = 0). 

Figureinta) shows the spectral function calculated with 
just the Fock component at different applied biases. We 
then increase the bias V by increasing the chemical 
potential of the left contact while keeping that of the 
right contact at zero. For low values oi V {V < uo) 
there is little change in the spectral function. However, 
once the value of V exceeds that of ujq, the spectral 
function becomes increasingly modified, especially when 
there is significant spectral weight inside the bias window 
fifi < uj < hl- In particular, the lineshapes of both the 
main peak and the vibron side-band peaks become de- 
formed, with a noticeable asymmetry of the main peak 
for biases where the main peak, but not the right-hand 
vibron side-band peak, lies within the bias window. 

There is a saturation regime once the vibron side band 
peaks, and thus nearly all of the spectral weight, is within 
the bias window (curves for /^^j ^1). Here the main peak 
becomes pinned around a; = 0.5 (effectively midway be- 
tween and /i/j), and a symmetric lineshape is restored. 
We postulate that this is owing to the bias being large 
enough to achieve simultaneous electron and hole trans- 
port. 

Figure [SJb) shows the spectral function with both the 




FIG. 6: Non-equilibrium spectral function for the off-resonant 
transport regime with applied bias for (a) the Fock diagram 
only and (b) both Fock and Hartree diagrams. The curves are 
offset vertically for clarity. The applied bias is given by the 
chemical potentials of the left (left-pointing arrows) and right 
(right-pointing arrows) leads respectively. The Hartree poten- 
tial has a strong effect on the peak positions. The other pa- 
rameters are eo ~ +0.5, 70 = 0.21, cjq — 0.3, toL,R = 0.15, 77 = 
0.005, »7y = 1. 



Hartree and Fock diagrams included. As we have already 
noted, at zero bias there is no change from the Fock-only 
spectral function. For small biases {V < ojq) there is 
little difference, but for larger biases the effect of the 
Fock diagram becomes increasingly evident. This is be- 
cause once the bias exceeds the vibron frequency wqj the 
non-equilibrium electron density becomes strongly per- 
turbed. The main peak is more stable in position than 
for the Fock-only spectral function, although it shifts 
slightly towards lower energies with increasing bias be- 
fore stabilising. The righthand vibron side-band peak 
becomes strongly deformed, and moves to lower, rather 
than higher, energies as the bias is increased. The left- 
hand vibron side-band peak appears at a lower bias than 
for the Fock-only spectral function, at a frequency just 
above zero, then tends towards zero as the bias increases. 
As for the Fock-only case however, once both vibron side- 
band peaks are within the bias window, saturation is 
reached and the peak positions stabilize. For both the 
Fock-only and Hartree-Fock spectral functions, the sepa- 
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FIG. 7: (color online) Hartree potential versus applied bias. 
The Hartree potential is almost constant in the low bias 
regime, and hence it could be neglected in the calculations, 
because it just corresponds to a global energy reference shift. 
However for larger biases, it is much more influential as the 
non-equilibrium electron density varies a lot with the applied 
bias. The parameters are the same as in figure [6l 



ration between both side-band peaks and the main peak 
is ~ 0.34, > Wo for the parameters we used. 

Figure [7] shows the Hartree potential (i.e. the real part 
of Sg_vib) plotted against applied bias for different values 
of the electronic level Sq. The solid line corresponds to 
the value of Eq used to calculate the spectral functions 
of the lower panel in Figure IHl The Hartree potential is 
small and almost constant in the low bias region. Hence, 
in the quasi-equilibrium regime, such a potential could 
be neglected in the calculations, because it simply corre- 
sponds to an energy reference shift. However when the 
bias window [/i^ , fij^] starts to encompass features in the 
spectral functions (either the main peak or vibron side- 
band peaks), the corresponding non-equilibrium electron 
density starts to vary substantially with the applied bias. 
Thus the Hartree potential shows a strong dependence on 
the value of the bias (as also shown in Ref. 60]), until 
the saturation regime is reached and the Hartree poten- 
tial is once again constant for very large biases. The 
asymptotic saturation value of the Hartree potential is 
dependent on the value of the electronic level Sq as one 
would expect. The first drop in the value of the Hartree 
potential happens at Sq, with the subsequent steps, which 
become progressively smaller and broader, at eo+ 0.34 
and then Eq + rtcjQ. 

This behaviour indicates that one cannot in principle 
neglect the Hartree diagram contribution in the calcu- 
lations, unless one is interested in calculating only the 
properties of the system for a bias range for which the 
electron density is (almost) constant^'^. 

It is worth mentioning that it is not straightforward 
to relate the modification of the peak positions in the 
two panels in figure [5] to the Hartree potential alone. 
In self-consistent calculations, a highly non-linear sys- 
tem needs to be solved, since the Hartree potential is ob- 




FIG. 8: Non-equilibrium spectral function A{lj) for the reso- 
nant transport regime for different applied biases for (a) the 
Fock diagram only, and (b) the Fock and Hartree diagrams. 
The curves are shifted vertically (-1-10 in y-axis) for clarity, 
with the values of the left and right chemical potentials given 
by the triangular symbols, here we have a symmetric potential 
drop. Adding the Hartree diagram breaks the electron-hole 
symmetry. The other parameters are the same as for Figure 
[S] except for eo = +0 and rjv — 0.5. 



tained from one element of the Green's functions which 
are themselves dependent on the values of the Hartree 
potential. 



2. Non-equilibrium resonant transport 

We now consider what happens when we apply a bias 
in the resonant transport regime. We apply a symmetric 
potential drop (i.e. the potential of the left electrode is 
raised by an amount rjeV while that of the right elec- 
trode drops by the same amount). This allows us to keep 
the electron-hole symmetry and see under which circum- 
stances the electron-hole symmetry is broken. 

The spectral functions for the Fock diagram are shown 
in figure [5Ja). As we have no Hartree term, the spec- 
tral function is electron-hole symmetric at equilibrium. 
Moreover, as we have chosen a symmetric potential drop 
(rj = 0.5) the spectral functions stay symmetric for all 
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applied biases. Constraining the symmetry in this way 
allows us to concentrate on the effects of increasing the 
bias, and we can see that the most important modifi- 
cations of the spectral function are in the width of the 
central and satellite peaks. The vibron side-band peaks 
show only a very small change with increasing applied 
bias. 

Figure ISI^b) shows the spectral function when we add 
in the Hartree term while keeping all other variables un- 
changed. We note immediately that the electron-hole 
symmetry is broken, even at equilibrium. In addition the 
main peak is shifted towards negative energies with re- 
spect to the Fock-only calculation, and the right-hand vi- 
bron side-band peak is completely suppressed — the whole 
spectral function takes on the qualitative appearance of 
a spectral function in the hole-dominated regime (com- 
pare with Figure U). On increasing the bias, the width 
and shape of the left-hand vibron side-band peaks is var- 
ied, and for higher biases new vibron side-band peaks 
appear above the main peak, until at high biases the 
spectral function is very similar in shape to that for the 
off-resonant spectral function (Figure IH]). 



C. A comparison with inelastic scattering 
techniques 

In this section, we will check the validity of the 
self-consistent Born approximation (SCBA), i.e. self- 
consistent calculations using only the Hartree and Fock 
diagrams, versus another method which gives more ex- 
act results as far as electron- vibron coupling is concerned. 
Since NEGF is a many-body perturbation expansion the- 
ory, by definition it does not contain all diagrams even 
though self-consistency allows us to achieve a partial re- 
summation of a subclass of diagrams. 

As an alternative to NEGF, one can calculate the 
transport properties using an extension of conventional 
scattering theory to include the interaction of incoming 
single-particle states with some bosonic degrees of free- 
dom inside the central region of interes t^^'^'^i'^^ . 

This technique, termed the multi-channel inelastic 
scattering technique (MCIST) ^^i^^ , has the advantage of 
being an exactly solvable problem, even in the presence of 
many non-interacting electronic states coupled to many 
vibration modes in the central region^^. MCIST is based 
on many-body perturbation theory for polaron and it is 
exactly solvable in the sense that MCIST is treating the 
electron-vibron coupling to all orders. In the language of 
polaron theory, MCIST contains all the diagrams corre- 
sponding to the electron-vibron interaction (in the corre- 
sponding transport regime). To be more precise, MCIST 
contains all orders of crossing and non-crossing diagrams 
in terms of the bare vibron propagato r^^i "^^^ 3^, In con- 
ventional polaron theory for one electron, there are no 
diagrams with electron-hole loops in them^l. 

However, MCIST is a single- particle scattering tech- 
nique and treats the statistics of the Fermi seas of the 



left and right leads only in an approximate manner. In 
the language of NEGF, this means that the results given 
by MCIST are only valid for a specific transport regime 
(as we will see below). 

In the case of the SSSM model, the retarded 
Green's function of the central region has the usual form 

G^iiu) = [ gliiu)-^ - - K-.M]-' • (28) 

Within MCIST, the retarded electron-vibron self-energy 
containing all orders of the electron-vibron coupling is 
expressed as a continued-fraction as shown analytically 
in Refs. ^ and 

ij 

Go(w-wo) 

GUu^ - 2a;o)-i - — 

G[;(w - 3wo) ^ - ■■• 

(29) 

where we recall that Gq is the retarded GF of the central 
region connected to the leads, i.e. Gq{uj) = [ 5q(cj)^^ — 
^rcads('^)]~^' without electron-vibron coupling. 

As show in Ref. [g^l, the lowest Born approximation, 
fully consistent with SCBA (Hartree- Fock calculations as 
shown above), is given in MCIST by an electron-vibron 
self-energy equivalent to Eq. (|29l) but where the integer 
n factors at each level of the continued fraction are all 
replaced by the integer n = 1. 

Obviously, this approximate substitution is good 
enough in the (very) weak electron-vibron coupling for 
which only the first level of the continued fraction con- 
tributes the most^. 

Below we compare the spectral functions obtained by 
exact MCIST calculations, by MCIST approximated to 
the SCBA (Hartree-Fock) level, and to NEGF-SCBA cal- 
culations. 

The spectral functions obtained from MCIST are 
shown in figure El for the weak/intermediate electron- 
vibron coupling regime (70/wo — 0.5) and for the strong 
electron-vibron coupling regime (70/wo = 0.8). The 
overall lineshapes correspond to a main peak with vi- 
bron side-band peaks located only above the main peak. 
These are typical results fully consistent with an off- 
resonant transport regime situation. Furthermore, the 
spectral functions obtained with BA-based approxima- 
tion (MCIST-SCBA and NEGF-SCBA) are virtually 
identical, especially in the intermediate (70 /wq = 0.5) 
to very weak (not shown here) electron-vibron coupling 
regime. For strong electron-vibron coupling, one obtains 
the same peak positions, however the amplitude of the 
peaks (especially the main peak) is slightly different. 

Since MCIST calculations will always give similar line- 
shapes independent of the value of £0, one can con- 
clude that MCIST calculations are only valid for the 
off-resonant transport regime at and near equilibrium. 
MCIST is not able to reproduce the spectral features 
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FIG. 9: (Color online) Equilibrium off-resonant spectral func- 
tions calculated with the MCIST technique and NEGF-SCBA 
for the off-resonant regime. The upper panel shows results 
for weak electron-vibron coupling ['yo/uio = 0.5), while the 
lower panel shows results for strong electron-vibron coupling 
(70/tJo = 0.8). MCIST calculations give exact results in this 
case, and can also be performed at the same level of approx- 
imation as SCBA (see main text for details). MCIST-SCBA 
and NEGF-SCBA are virtually identical, but show large dis- 
crepancies for the peak positions (main peak and more espe- 
cially for the vibron side-band peaks) in comparison to the 
exact results. The parameters are eo = -1-0.5 (electron trans- 
port), 70 = 0.2 (upper panel) and 70 = 0.32 (lower panel), 
a;o = QA,toL,R = 0.15. 



of the resonant transport regime (i.e. vibron side-band 
peaks on both sides of the main peak). This is essentially 
due to the fact that in MCIST, one does not take prop- 
erly into account the statistics of the Fermi seas of the 
left and right leads. See for example Eqs. (^5)) and ([^. 
there are no leads' Fermi distributions in the retarded 
component of the leads' self-energies I]£ ^j. 

Now, comparing BA-based calculations with the exact 
MCIST calculations, one can see from figure |9] that the 
BA-based calculations give the wrong polaron shift, i.e. 
the normalized position of the main peak Eqj especially 
in the strong electron-vibron coupling regime. Further- 
more BA-based calculations also give the wrong energy 
separation between the main peak and the first vibron 



side-band peak. This energy difference should be equal 
to the vibron energy wq, as it is given by exact MCIST 
calculations. Note that the limits of BA-based calcula- 
tions were also been studied by Lee et al. in a somewhat 
different context in Ref. 60]. 

In conclusion, this means that Hartree-Fock (or BA) 
based calculations for electron-vibron interaction are only 
valid for weak coupling, as can be expected from a 
perturbation-expansion based theory. Hence one needs to 
include higher-order diagrams in the electron-vibron self- 
energies to go beyond the commonly-used self-consistent 
Born approximation (Hartree-Fock) in order to obtain 
correct results for a wide range of parameters. The ef- 
fects of the higher-order diagrams (here second-order- DX 
and DPH diagrams) are explored in detail in the follow- 
ing sections. 

Additionally, although MCIST calculations are only 
valid in the off-resonant transport regime at and near 
equilibrium, they include all possible higher-order dia- 
grams (with bare vibron propagator) and hence can be 
used as a reference for any perturbation-expansion-based 
NEGF calculations performed at equilibrium or in the 
quasi-equilibrium regime. 



D. Vertex corrections and polarization effects to 
the spectral functions 

In this section, we present results for the spectral func- 
tions when the second-order diagrams (see Figure [3]) are 
included in the calculations of the Green's functions. The 
reader can find more information about the mathemati- 
cal expressions for the self-energies corresponding to the 
second-order diagrams in Appendix 1X1 

These diagrams fall into two types — the double- 
exchange DX diagram, corresponding to vertex correc- 
tions, and the dressed vibron diagram, which includes 
a single electron-hole bubble, renormalizing the vibron 
propagator and hence giving rise to polarization effects. 

We have used three different levels of approximation 
to calculate these Green's functions: Firstly, calcula- 
tions with no self-consistency — the Green's functions are 
simply calculated using the diagrams in Figure [5] and 
Figure [3] using the bare propagator Go as the electron 
Green's function. In our model, Gq is the Green's func- 
tion of the central region connected to the leads with no 
electron-vibron interactions. This is a first-order pertur- 

, ■ r u ■ u ^H,F.DX, and/or DPr^ i 

bation expansion for which ' [CJo\. We 

use the abbreviations BA (Born approximation for non 
self-consistent Hartree and Fock diagrams) and BA-I-DX 
(DX for double exchange) and BA+DX+DPH (DPH for 
dressed vibron, for the GMK-like diagram) in the follow- 
ing. 

Secondly, we perform partly self-consistent calcula- 
tions, where the Green's functions are calculated with the 
first loop of self-consistent calculations with the Hartree 
and Fock diagrams. We use these Green's functions as a 
starting point to calculate new, corrected. Green's func- 
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tions including the second-order diagrams S^^^j,^ and/or 
^?-^h [G^'~^^^] in which the Green's functions are the cor- 
responding SCBA Green's functions. 

Finally, we perform fully self-consistent calculations, 
in which the Green's functions are calculated in a 
self-consistent manner with the first and second-order 
diagrams included within each iteration of the self- 
consistency loop. 

Our rationale for choosing to do calculations in this 
manner is as follows. We are able to test the different lev- 
els of approximation and see very precisely the effects of 
vertex corrections and polarization are on both the bare 
Green's functions Go and the SCBA-level Green's func- 
tions. We will also show in section UlI D II that by using 
the SGBA Green's functions rather than Go as a start- 
ing point, one achieve a better convergence in the cal- 
culations. In particular we will show later on in Section 
IIIID II that in some cases, for example the off-resonant 
transport regime when 70 /wo > 0.6 — 0.7), the use of the 
bare Go Green's functions as a starting point for a fully 
self-consistent calculation with second-order diagrams ac- 
tually gives unphysical results. Additionally, fully self- 
consistent calculations are extremely computationally in- 
tensive, and hence it is both interesting and useful to ex- 
plore the range of parameters for which the second-order- 
diagram corrections to the SGBA calculations give a suf- 
ficiently accurate description of the Green's functions in 
comparison to the fully self-consistent calculations. 



1. Off-resonant regime at equilibrium 

Figure [TO] shows the spectral functions of the off- 
resonant transport regime at equilibrium, calculated for 
weak/intermediate electron- vibron coupling (70/wo = 
0.50) and for different diagrams and levels of self- 
consistency. 

As already mentioned in sections IIII Al and IIII Bl the 
energy separation between the main peak and the first 
vibron side-band peak (which should be equal to the vi- 
bron energy, here wo = 0.40) is not well reproduced by 
BA-based (Hartree-Fock based) calculations and is much 
larger than ujq. The self-consistency introduced in the 
calculations give a marginally smaller energy separation. 

The effects of the second-order DX diagram are firstly 
to bring the vibron side-band peak closer to the main 
peak, hence giving an energy separation closer to the ex- 
act value; and secondly, one observes a strong narrowing 
and a larger amplitude of the vibron side-band peak in 
both self-consistent and non-self-consistent calculations. 
Both these effects thus qualitatively modify the spectral 
functions towards better agreement with the exact results 
as shown in Section [III CI 

It is, however, worth mentioning that including the DX 
diagram in the calculation does not greatly affect the po- 
sition of the main peak, a result which may be under- 
stood from the fact that the dynamical polaron shift is a 
quantity difficult to obtain exactly from a perturbation 
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FIG. 10: (Color online) Equilibrium spectral function 
for the off-resonant electron transport regime in the 
weak/intermediate electron- vibron coupling. Non self- 
consistent (top panel) and self-consistent (bottom panel) 
calculations were performed with first-order Hartree and 
Fock-like (BA, solid lines), plus second-order DX (dashed 
lines), plus DPH (squares) electron- vibron self-energy dia- 
grams. The other parameters are eo ~ +0.5, 70 ~ 0.2, cjo = 
0.4:,toL,R = 0.15, — 0.025. Self-consistent calculations give 
better spectral functions, with vibron side-band peak closer 
to the main peak. The second-order DX diagram narrows the 
side-band peaks and bring them even closer to the main peak 
as expected. The second-order DPH diagram does not con- 
tribute much in the off-resonant regime. In this regime and 
for small electron-vibron coupling, second-order corrections to 
SCBA calculations are a good approximation to correspond- 
ing full self-consistent calculations. 



expansion theory beyond weak coupling^iS. 

The effects of the second-order DPH diagram in 
the spectral function are virtually nil for the case of 
weak/intermediate coupling and the off-resonant regime. 
This might not be that surprising since the DPH diagram 
corresponds to a Fock-like diagram with a renormalised 
vibron propagator. The renormalization of the vibron 
is due to a single electron-hole bubble. However in the 
off-resonant transport regime, the spectral function is al- 
most empty (in the case of electron transport) or almost 
full (in the case of hole transport) which implies that 
there are not many electron-hole excitations available in 
this transport regime. Hence the polarization (i.e. the 
contribution of the electron-hole bubble) is very small, 
subsequently giving very small values for the DPH self- 
energy. As an example, we have checked our numerical 
values in the case shown in figure [TO] and we have found. 



as expected, that the maximum values of E 



DPH,i 
-vib 



are 30 



DX,r 
e-vib ' 



to 50 times smaller than the maximum values of S 

Furthermore, in the off-resonant regime and for 
weak-ish electron-vibron coupling, it seems that the 
fully self-consistent (full SC (BA-fDX) curve in fig- 
ure nop results Gj^jj gc(w) are well approximated 
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by the results given by a second-order correction 
to a self-consistent Hartree-Fock (SCBA) calculations 
(SCBA+DX curve in figure m G;,,,,,, sci^) = 

[G'scba('^)"^ - ^^^■'■[Gscba] ]~^- This is an interest- 
ing result as it implies that the physical properties of the 
system, at least in these conditions, could be well de- 
scribed by a second-order correction of the lowest-order 
SC calculations, without the need to perform a fully self- 
consistent calculation up to the second-order. 

Now we turn to the analysis of results obtained for 
stronger electron-vibron coupling. Examples of such cal- 
culations are given in figure [Til for 70/1^0 = 0.80. In the 
strong electron-vibron coupling regime, one obtains qual- 
itatively the same contributions of the second-order dia- 
grams as explained above for the weak coupling regime. 
The DPH diagram does not have a large role in the off- 
resonant transport regime, although slightly affecting the 
width of the main peak. And the DX diagram shifts the 
vibron side-band peak towards the main peak (and hence 
towards the exact results) as well as narrowing the peak 
width and increasing the peak amplitude. These effects 
are amplified in Figure [11] because the electron-vibron 
coupling constant 70 is bigger than in Figure [TUl Fur- 
thermore, for strong electron-vibron coupling, the DX 
diagram modifies the energy position of the main peak 
and seems to give a slightly better polaron shift. 

There is however another interesting effect observed 
from such a set of calculations, which can be seen in the 
upper panel of Figure [TT] The first iteration of a fully 
self-consistent calculation (including the Hartree, Fock, 
and DX diagrams) starting with the non-interacting GFs 
Go gives non physical results, i.e. negative values of 
the spectral function (see curve B A-I-DX in figure [TlT) , 
and such an unphysical behaviour does not self-correct 
in the following iterations. This a well known problem, 
which has already been encountered in the past by sev- 
eral authors in the context of electron-vibron interac- 
tion (for example in Ref [t^) and also in the context 
of electron-electron interaction when considering topo- 
logically equivalent diagramsSiSflSl. In Ref. [t^ , a some- 
what different approach than ours was used. It is based 
on a linked cluster expansion for the non-equilibrium 
steady-state regime, and negative densities of states were 
obtained when including second-order diagrams, and in 
some cases even higher-order cluster approximations did 
not seem to give a convergent solution at intermediate 
electron-vibron strength and in the presence of the Fermi 
seas. 

However, our calculations reveal that it is possible to 
solve such a problem by starting the fully self-consistent 
calculations (up to second order) from a different starting 
point, namely by starting from the Gscba Green's func- 
tions (i.e. the GFs obtained from a fully SC calculation 
including only the lowest-order diagrams) . This is shown 
by the curve (SCBA-I-DX) in the upper panel of figure [TT] 
and by the lower panel in which all the fully converged 
self-consistent results are shown. 

For the moment, we do not have a full physical expla- 
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FIG. 11: Equilibrium spectral function for the off-resonant 
electron transport regime and strong electron-vibron cou- 
pling. Calculations were performed with different dia- 
grams: first-order Hartree and Fock-like (BA) and second- 
order (DX, DPH) electron-vibron self-energy diagrams. (Top 
panel): Calculations were done with no self-consistency 
(BA,BA+DX), partial (SCBA+DX) and full self-consistency 
(SCBA). (Bottom panel): full self-consistent SC(BA+...) cal- 
culations. The other parameters are eo ~ -1-0.5, 70 = 
0.32, a;o = 0.4, foi,^ = 0.15, »7 = 0.025. The calculations show 
that no self-consistency gives poor results in comparison to ex- 
act calculations given in Figure [5] A closer spectral function 
to the exact result is obtained from SC(BA-I-DX) calculation. 
The second-order DPH diagram do not contribute much in 
the off-resonant regime. See main text for a detailed analysis 
of the spectral functions. 



nation of the reason why starting from a SCBA calcula- 
tions is better to achieve full self-consistency with higher- 
order diagram than Hartree-Fock, apart from the simple 
fact that an Hartree-Fock calculation is probably closer 
to the true interacting solution than the non-interacting 
solution. 

To conclude this section, we can say that in the off- 
resonant regime at equilibrium, the second-order DX di- 
agram dominates over the second order DPH diagram. 
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FIG. 12: (color online). Equilibrium spectral function for the 
resonant transport regime and weak/intermediate electron- 
vibron coupling. Calculations were performed with diagrams: 
BA (Hartree-Fock) and second-order DX, DPH electron- 
vibron diagrams. The other parameters are eo = +0.0 (res- 
onant transport), 70 = 0.24, tJo = 0.4, toi.fl = 0.2, 77 — 0.03. 
The top inset shows the whole spectral functions, while the 
left (a) and right (b) parts are zooms, for positive u, of the 
central peak and first vibron side-band peak respectively. The 
second-order diagrams have the following effects on the SCBA 
spectral functions: DPH broadens the central and side-band 
peaks as well as lowers their height, in opposition to DX which 
narrows the peaks and brings the side-band peaks slightly 
closer to the center, with a strong and unphysical increase of 
the amplitude of the central peak. Calculations with both 
DX and DPH give in appearance a broadening of the central 
peak, which however recovers the correct height necessary to 
conserve the Fermi-liquid property of the electron-hole sym- 
metric system (see Section rill El for more detail). 



2. Resonant regime at equilibrium 

In this section, we present calculations for the resonant 
regime at equilibrium. Though we have shown, in sec- 
tions nil A 21 and HUB 21 the importance of the Hartree 
diagram, we will consider below results obtained with- 
out the Hartree diagram. The calculations were per- 
formed with the Fock and second-order DX and/or DPH 
diagrams which conserve (at numerical accuracy) the 
electron-hole symmetry of the system. 

The reasons why we have choosen to perform this 
model calculation are twofold: first, even for a single 
electronic level, we expect to have the maximum possi- 
ble electron-hole excitations available when the spectral 
functions are electron-hole symmetric. Hence we expect 
the polarization effects be to more pronounced in a sys- 
tem with electron-hole symmetry. Second, the electron- 
hole symmetric model permits us to emphasize the com- 
petitive effects between the DX and DPH diagrams as 
will be shown below. 

Figure [T^ shows the equilibrium spectral functions 
of the resonant regime (at weak/intermediate electron- 
vibron coupling) obtained from different self-consistent 
calculations including first and second-order diagrams. 



On one hand, the second-order DPH diagram corre- 
sponds to a partial dressing of the vibron propagator 
by one electron-hole bubble. It gives an extra lifetime 
in the retarded/advanced electron Green's functions in 
comparison to Hartree-Fock (SCBA) calculations. Hence 
the main effect of the DPH diagram is to introduce an 
extra broadening of the peaks in the spectral functions. 
Since the DPH diagram is one of the so-called conserving 



approximation 



.71.81 



the broadening of the peaks leads 



to a reduction of their height to keep globally the same 
total spectral weight. These effects can be seen on the 
SC(BA+DPH) curve in Figure [H 

On the other hand, the second-order DX diagram, 
which is a conserving approximation, has an opposite 
effect: a strong narrowing of the peaks (especially of 
the central peak) accompanied with an increase of their 
height, as can seen on the SC(BA+DX) curve in Figure 

Now it is interesting to see what happens when 
the calculations are performed with both second-order 
DX and DPH diagrams. This is shown on the 
SC(BA-|-DX+DPH) curve in Figure [H We obtain a hy- 
brid behaviour, in the sense that in appearance the cen- 
tral peak is broadened in comparison to SCBA calcula- 
tions. However, the height of the peak at the Fermi level 
is conserved (up to numerical accuracy). This is a very 
important result which proves that for the model calcu- 
lation of an electron-hole symmetric system, one has to 
include both the DX and DPH second-order diagrams in 
order to conserve the expected Fermi-liquid properties of 
the system. In this regime both the DX and DPH second- 
order diagrams play an equally important role which de- 
termines the linear response properties of the system as 
shown below in Section UlI El 

By comparison, in a diagrammatic treatment of the 
electron- electron interaction on the electron propagator, 
the situation is often different: the electron- hole bub- 
ble diagram that appears here in the DPH contribution 
to the propagator is large, especially in highly polariz- 
able metallic and open-shell systems where electron-hole 
pairs may be created with low energy cost, because the 
Coulomb interaction operates at all energy scales. In that 
case, summing the bubble diagrams to infinite order as 
is done in Hedin's GW approximatior^iS- is much more 
important than including the second-order exchange dia- 
gram. The key difference for the electron-phonon interac- 
tion is that the vibron frequency luq imposes a restricted 
energy scale on the interaction, reducing the importance 
of the bubble diagrams, and correspondingly increasing 
the importance of the second-order exchange diagram. 



3. Off-resonant regime at finite bias 

Figure [T3] shows the spectral functions of the off- 
resonant transport regime at finite bias. The calculations 
have been performed for weak/intermediate electron- 
vibron coupling (jq/luq = 0.50) and for different dia- 
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grams and levels of self-consistency. We have considered 
the case for which real excitations of vibrons are possible 
V = 0.45 > loq. 

The interpretation of the results is not as straightfor- 
ward as in the equilibrium case, because non equilib- 
rium effects are sometimes counter intuitive. However, 
the overall shapes of spectral functions are quite similar 
to that obtained within Hartree-Fock based calculations 
(see Section ITlIB 1|) in the sense that they present a cen- 
tral peak with vibron side-band peaks on both sides. 

Similarly to the resonant regime at equilibrium, the 
main effects of the second-order DX diagram is to narrow 
the width of all peaks, and to shift slightly the side-band 
peaks towards the main peak. 

As mentioned in the previous section, the main effect of 
the second-order DPH diagram is to broaden the peaks, 
in opposition to the effects of the DX diagram. How- 
ever here, the broadening in the non-equilibrium condi- 
tion appears less important than in the case of the res- 
onant regime at equilibrium. Hence full self-consistent 
calculations performed with both DX and DPH second- 
order diagrams give a narrowing of the peaks with a cor- 
responding increase of their amplitude in comparison to 
SCBA calculations. 

Finally, out of equilibrium, it can be seen that the 
results given by a full self-consistent calculation (curve 
SC(BA-I-DX) in Figure [T3)) are strongly different from 
a second-order correction to an Hartree-Fock calculation 
(SCBA-I-DX curve in[T^ which gives an excessively nar- 
rowed central peak. 

This means that in the weak/intermediate electron- 
vibron coupling, second-order corrections to a SCBA cal- 
culations are only good enough at equilibrium, however 
at non-equilibrium full self-consistency needs to be per- 
formed with all diagrams of the same order. 




(0 

FIG. 13: (color online). Non-equilibrium spectral function 
(bias V — 0.45 > too) for the off-resonant transport regime 
for the weak/intermediate electron-vibron coupling. Calcula- 
tions were performed self-consistently with first-order Hartree 
and Fock-like diagrams (SCBA curve), and with second-order 
DX and/or DPH diagrams (SC(BA+DX,+DPH,+DX+DPH) 
curves). Second-order DX correction to SCBA is also shown 
(SCBA-I-DX curve). The top- left inset shows a zoom on 
the central peak around cj ~ 0.4, and the top-right in- 
sets show a zoom on the vibron side-band peaks around 
uj ~ 0.0 and uj ~ 0.8. The other parameters are £o = -1-0.5 
, 70 = 0.2, ojo = OA,toL,R = 0.15,7? = 0.025, rjv = 1. The 
second-order DX diagram narrows the peaks with a slight 
shift of the side-band peaks towards the central peak. The 
narrowing is too strong in the case of partially self-consistent 
calculations (SCBA-I-DX). The second-order DPH diagram 
broadens the peak, but not as much as in the resonant case. 
Full self-consistent calculations including both second-order 
diagrams result in an intermediate behaviour for the modifi- 
cations of the spectral functions. 



4- Resonant regime at finite bias 

As we have already shown in the first-order electron- 
vibron diagrammatic calculations, the spectral func- 
tions for the off-resonant and resonant regime at non- 
equilibrium are qualitatively similar, in the sense that 
they present a central peak with vibron side-band peaks 
on both sides. One can compare for example the spec- 
tral functions obtained for Hartree-Fock-like calculations 
at non-equilibrium shown in Figure [5] and Figure [S) 

Hence, and again on a qualitative level, the effects of 
the second-order DX and DPH diagrams on the resonant 
case at finite bias are similar to what has been obtained 
for the off-resonant non-equilibrium case described in the 
previous section. The effects of these higher-order dia- 
grams on the full non-equilibrium transport properties 
for the different transport regimes will be presented in a 
forthcoming paper. 

However before turning the discussion to the linear- 
response properties of the system at and near equilib- 
rium, we would like to comment on a specific aspect of 



the effects of higher-order diagrams. 

The narrowing of the vibron side-band peaks and of the 
main central peak due to higher-order (DX) diagrams was 
also obtained by other authors (see for example Ref. 84] ) . 
In this paper, a different non-equilibrium approach was 
used. It consists of starting with an electron dressed by 
a vibron (a polaron) in the isolated central region, then 
using perturbation expansion theory (with partial resum- 
mation) in terms of the coupling of the central region to 
the non-equilibrium left and right leads. However, the 
results for the spectral functions in Ref. [Sj were only 
given for the resonant transport regime at equilibrium 
and all calculations were performed without the Hartree- 
like diagram. Though we have already shown that such 
a diagram plays a crucial role in the spectral properties 
of the system in the resonant regime at and out of equi- 
librium. 



15 



tjlin/(j-0 


off-resonant 


off-resonant 


resonant 




weak e-vib 


strong e-vib 


(e-h symmetric) 




r Ig. IIUI 


i'lg. [iij 


i'lg. [iZJ 


SCBA 


0.002654 


0.004470 


0.8379* 


bOBA+Da 


0.002666 


r\ r\r\ a r'r\r' 

0.004606 




SCRA+DX+DPH 


OOSlSfi 






SC(BA+DX) 


0.002668 


0.004606 


1.1296 


SC(BA+DPH) 






0.7194 


SC(BA+DX+DPH) 




0.004638 


0.8545 


MCIST 


0.002621 


0.004936 


n.a. 


no e-vib 


0.002021 


0.002021 


1. 



* in principle here G\\^ should be Go but is not because we are 
using a tiny but finite r\ = 0.030 value. See appendix IbI for 
detailed explanations. 

TABLE I: Linear conductance Gun = dl /dV\v-^o obtained 
from the value of the equihbrium spectral density at a; = fi""^. 
Giin is given in units of the quantum of conductance Go ~ 
2e^//i. The equilibrium spectral functions are shown in Figs. 
1101 111! and I12p and were done for different different lev- 
els of approximation as explained in the corresponding figure 
captions. 



E. Linear response transport properties 

We now briefly discuss the effects of electron-vibron 
interaction taken at different levels of the diagrammatic 
expansion on the linear-response transport properties of 
the single-molecule nanojunction. Before doing so, we 
explain how to derive the linear conductance from the 
value of the spectral functions calculated at equilibrium. 

In our model, there is a direct proportionality between 
the left and right leads' self-energies S£ ^ because, the 
central region is coupled to the leads via the single hop- 
ping matrix element and we have chosen identical leads. 
Then the current (Eq. ([^0]) ) can be recast as follows (for 
the details of the derivation, see for example Refs. |64ll7g |): 

/ = -4e/h f duj ihioj) - fnioj)) r{^)^m[G^ {uj)], (30) 



where T = TLTR/iTL + Tr) and 
The linear conductance 



-25m[S^] 



lin 



dl 

dV 



is obtained from 



Gn„/Go = -2T{ti'^^) 3m[G'-(Ai"'i)] = 2r(M°q) Aifi^'^) 

(31) 

where Gq — jh is the quantum of conductance. 

Table IIII El shows the different values for the linear 
conductance Gnn obtained from the equilibrium spectral 
functions shown in figures [TOl [TT] and [121 

For the off-resonant regime at the weak electron-vibron 
coupling, there are not many differences in the linear 
conductance Gun values calculated for all the different 
diagrams, as one might expect from the spectral func- 
tion behaviour shown in figure [TUl This is especially true 



when the main peak in the spectral function is located 
well away from the Fermi level, i.e. |eo~M°''l ^ linewidth 
of the peaks. For strong electron-vibron coupling, the 
differences between different levels of approximations are 
more pronounced. All Gun values are smaller than the 
exact MCIST result. This is essentially due to the fact 
that perturbation expansion gives only an approximate 
value of the polaron shift. The differences from the exact 
result are more important for SCBA-based calculations 
(difference of ~ 10% with exact result) and of course 
calculations performed up to second order give a better 
polaron shift and hence better values for the linear con- 
ductance (~ 5% difference only). 

In the quasi- resonant case when jeo — '^-^ (1 to 
2) linewidth, the effects described above will be more 
pronounced, because the linear conductance is no longer 
given by the extreme tail of the main peak crossing the 
Fermi level. 

For the resonant case in the absence of the Hartree po- 
tential, the SCBA-based calculations should give a per- 
fect linear conductance. The reason why Gun 7^ Gq is 
because the calculations were performed with a tiny but 
finite value of r] as explained in detail in Appendix |B] 
Hence one may say that, for the chosen set of parame- 
ters, the value Gun = 0.838 represents an upper bound 
for the linear conductance (the corresponding numerical 
perfect conductance). And as explained in Appendix iBl 
linear conductance values can only be compared between 
calculations performed with the same value of 77, which 
is what we have done. 

The linear conductance is strongly renormalized (de- 
creased and increased) when including only one of the two 
second-order diagrams (DPH and DX respectively). In 
contrast, it becomes close again to the expected quantum 
of conductance when the calculations are performed with 
both second-order diagrams. This dependence of the con- 
ductance is well understood from the behaviour of the 
spectral functions shown in Figure [T2l SC(BA-I-DPH) 
calculations introduce an extra broadening and a cor- 
responding decrease of amplitude of the peaks in the 
spectral functions, hence a decrease of Gun- However, 
the strong narrowing, with increased amplitude, of the 
peaks due to SC(BA+DX) calculations leads to an in- 
crease of Giin. The SC(BA+DX) calculations give a lin- 
ear conductance larger than that obtained from SCBA, 
i.e. Giin > Go, which is an unphysical result for our 
non-degenerate electronic level in the central part of the 
system. 

This implies a strong constraint on the validity of the 
values of Gun obtained from a many-body perturbation 
expansion: in order to conserve the Fermi-liquid prop- 
erties (see Appendix |B]) in the electron-hole symmetric 
resonant regime, one has to perform the calculations by 
including all diagrams order by order. Results obtained 
from partial resummation of a subset of diagrams will 
probably give an incorrect linear conductance. 

Hence calculations performed by fully renormalizing 
the vibron propagators^M^ or all crossing diagrams^. 
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while containing diagrams higher than second order, will 
probably break important physical properties of the sys- 
tem at T = 0, because they lack important electronic 
processes, to second-order the DX vertex correction-like 
diagram or DPH diagram respectivelly. 

In the resonant (or quasi-resonant) regime, polarisa- 
tion effects and vertex corrections play an equally impor- 
tant role in the electronic structure and transport prop- 
erties of electron- vibron interacting nanojunctions which 
electron-hole symmetry. 



IV. CONCLUSION 

In this paper, we have presented a method based on 
NEGF to calculate the equilibrium and non-equilibrium 
electronic structure of electron- vibron interacting single- 
molecule junctions. We have applied the method to a 
model system which consists of a single electronic level 
coupled to a single vibration mode in the central region, 
the latter in contact with two non-equilibrium electron 
reservoirs. 

In comparison to previous studies performed within 
a similar approacb22*25^>28^-34^d4-46^ ^-^^ j^o^gl^y of 

our method lies in the fact that it goes beyond the 
conventionally-used self-consistent Born approximation 
(SCBA). Higher-order diagrams for the electron- vibron 
interacting have been implemented in our calculations. 

In this paper we have considered the second-order dia- 
grams which contain two-vibron processes: the so-called 
double exchange DX diagram which is part of the vertex 
corrections to SCBA, and the partially dressed DPH dia- 
gram in which the vibron propagator is renormalized by 
one electron- hole bubble. 

We have studied the effects of the first and second- 
order diagrams on the spectral functions for a large set 
of parameters and for different transport regimes (reso- 
nant and off-resonant cases) at equilibrium and in the 
presence of a finite applied bias driving the system out 
of equilibrium. 

We have shown the important role played by the 
Hartree diagram in calculations based only on first-order 
diagrams. Such a diagram should not be neglected un- 
less one works within an applied bias range for which the 
correspond Hartree potential is constant. 

For calculations including both first and second-order 
diagram we have found that the effects of the individ- 
ual second-order diagram are as follows: the DX dia- 
gram reduces the width of the peaks in the spectral func- 
tions, with an increase of their height, while the DPH has 
the opposite effect: it increases the width and decreases 
the height. Furthermore, the DX diagram moves the vi- 
bron side-band peak position towards the main peak, and 
hence gives a better peak separation as should be ob- 
tained from exact calculations. Calculations performed 
with both DX and DPH diagram give intermediate re- 
sults for the spectral functions, which are not simply an 
average/superposition of the individual effects because of 



the strong non-linearity involved in solving the problem 
self-consistently. 

Furthermore, the effects of the second-order diagrams 
also depend on the transport regime: at and near equilib- 
rium, the DX diagram dominates over the DPH diagram 
in the off-resonant transport regime (essentially because 
there are not many electron-hole excitations available in 
this regime). In the resonant case, however, both DX 
and DPH play an equally important role. For large non- 
equilibrium conditions {V > lu), both second-order di- 
agrams play an important role, since the corresponding 
spectral functions look qualitatively similar to those ob- 
tained the resonant regime. 

We have thus shown that it is indeed necessary to 
go beyond SCBA to obtain the correct results for a 
wide range of parameters. This has also been confirmed 
by comparing our NEGF results to an exact calcula- 
tion in terms of electron-vibron interacting (the multi- 
channel inelastic scattering technique MCIST^*), though 
it should be noted that comparison between NEGF and 
MCIST calculations are only valid for the off-resonant 
transport regime. 

We have also studied in detail the effects of self- 
consistency on the calculations, especially for the cases 
including the second-order diagrams (see Section fill D II 
and Figure [TT]). We have found a solution to an old prob- 
lem, well known in many-body perturbation theory: in 
order to avoid negative spectral densities when including 
higher order diagrams in the calculations, it seems more 
appropriate to start the self-consistent loop with SCBA 
(Hartree-Fock like) electron Green's functions . 

Finally, we have studied the linear response (linear con- 
ductance Giin) of the nanojunctions, and have found that 
in the off-resonant regime, the value of Gun is governed 
by the behaviour of the tail of the main peak at the Fermi 
level, hence it depends on both the position and the width 
of this peak. In the off-resonant regime, the DPH dia- 
gram contribution is negligible and the DX diagram gives 
a better polaron shift (better position of the main peak) 
and hence the second-order (mostly DX) calculation gives 
a better agreement for Gun with the exact result. For the 
near-resonant regime, the contribution of DPH will be- 
come more important and both DX and DPH will play 
an important role in determining the value of Gun by 
changing both the position and the width of the main 
peak. 

For the resonant regime with electron-hole symmetry, 
it is necessary to perform the calculations by includ- 
ing both DX and DPH diagrams in order to conserve 
the Fermi-liquid properties of the system. We anticipate 
that this will be also true for higher-order diagrams and 
calculations need to be performed by including all dia- 
grams of the same order. Calculations performed with 
partial resummation of a subset of diagrams will break 
the expected Fermi-liquid properties of the system and 
will probably lead to an incorrect value of the linear con- 
ductance. 

Finally, we expect that the effects of the second-order 
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diagrams, shown in this paper for the spectral functions 
only, are also important in the full non-equilibrium trans- 
port properties of the nano junctions. The study of these 
effects on the non-linear conductance are currently un- 
dertaken and will be considered in a forthcoming paper. 



Appendix A: Self-energies for the Born 
approximation and the second-order diagrams 

In this appendix, we show how to derive the expression 
for the self-energies corresponding to one- and two-vibron 
process diagrams in terms of electron-vibron coupling. 

Starting from the SSSM Hamiltonian, the general def- 
inition of the one-particle GF of the central region is ob- 
tained from the time- loop contour-ordered product of the 
creation and annihilation operators in the central region 

G(r,r') = -i(Tc,dff(r)4(r')). (Al) 

The time-loop contour Ck contains two branches, the up- 
per and the lower (— ) branch. On the upper branch, 
time starts in the infinitely remote past and evolves for- 
wards, then at the turning point, which can be placed 
at any arbitrary time, one passes onto the lower branch 
where the system evolves backwards in time back to the 
initially non-interacting starting point at t — —oo. 

Then any expectation value of products of opera- 
tor reduces to (0o|7ca' ■ ■ ■ ^Ck^ where 
((/)o| . . . |'/>o} is the average over the non- interacting 
ground state. The operators are then given in the 
interaction picture, and Sck is the generalization of 
the time evolution operator on the Keldysh contour 

Sck = Tc,, (exp{-i dr V{t)}^ where T^, is the 

time-ordering operator on the contour Ck and /^^ dt 

implies integration over Ck- And V is the "perturba- 
tion" to the reference Hamiltonian, which in our case 
would be the interaction electron-vibron as well as the 
coupling of the central region to the leads). 

Expanding Sck as a series in terms of the electron- 
vibron coupling Hamiltonian i?e-vib, one can derive the 
electron-vibron self-energies to any order of the electron- 
vibron coupling by calculating any time ordered products 
in the series using the usual rules of many-body perturba- 
tion theory, such as Feynmann diagrammatic expansion 
or Wick's theorem^^"^. 



1. Lowest order self-energies 

The one-vibron process self-energies, corresponding to 
the Hartree and Fock-like diagrams (see Figure ^ are 
given by 

Sffvib (n , r2 ) = -i7o' ^0 (ri , ) G(ri , t+ ) , ( A2) 



and 

5^f-vib (n 1 ) = i7o Do{ti,T2) G{ti,T2) , (A3) 

where are times on the time-loop contour. 

The projections onto the real (physical) times are given 
by, for example, 

^!^:^^^'{tut2) = 111 Di^^'ihM) G^^^^ihM) , (A4) 

for the Fock self-energy . 

The index Q = ± labels the branch of the time-loop 
contour corresponding to the forward (^i = +)/ back- 
ward (Ci = — ) time evolution respectively. 

There are several useful relationships between the 
different projections (or Keldysh components) of the 
Green's functions X = G, They are 

X"" ^ X++ - X+- = X-+ - X— (A5) 

or equivalently in terms of time-ordered {t = H — h), anti 

time-ordered (f = ), greater (>= — h) and lesser (<= 

H — ) Green's functions: 

X"^ = X* - X< = X> - X* (A6) 

and 

X" ^ X++ - X-+ = X+- - X— (A7) 

(or equivalently X"" ^ X* - X> = X< - X*), and hence 
-X" = X> - X<. 
Similar relationships exist for the self-energies O 

{x = r,a,<,>,t,i). 
Using these relationships and taking the steady state 
limit, i.e. X{t,t') = X{t — t'), and after Fourier trans- 
formation into an energy representation ^(a;), we ob- 
tain the usual expressions20*24^^^-^^i4-46 

Hartree T,f_':^-^^{uj) and Fock S^'^jj^(w) electron-vibron self- 
energies {x = r, a, <, >). 

For example, the Hartree and Fock self-energies are 

- -hi = 0) / ^ ^^('^'^ ' 
^e-Li^) = nl J ^ m^') G'i^ . (A9) 

where y represents one of the tree Keldysh compoentns 
y = i, <,>. Using the relationship E*" = — we 
find 

sr.vibH-i7o'/ ^i?5(---')G<(..') 

+ D<{oj~ tj')G''(t^') + DJ^ioj - ^^^(a;) . 

(AlO) 
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With the usual definitions for the bare vibron Green's 
functions Dq: 



D>{u;) 



+ ((iVph) + l)(5(c^ + c^o)] 
-2TTi [(iVph)(5(w + Wo) 

+ ((iVph) + l)<5(c^-wo)] 
1 



Ul — UlQ + IT] 
1 



UJ + UJo+ IT] 



(All) 
(A12) 

(A13) 
(A14) 



where (iVph) is the averaged number of excitations in 
the vibration mode of frequency ujq given by the Bose- 
Einstein distribution at temperature Tvib- 

One can see that T,fj^^^'^{uj) — because Dq'^{uj — 
0) = unless wq = which would be an odd case of 
study. And since Do'°(w = 0) = — 2/cjo, one has 

sf;ib = sf;tb-^2i|^G<(.'). (A15) 

Furthermore the lesser and greater Fock self-energies 
can be expressed in a more compact form: 



"'e-vib 



M=7o [(^ph)G^<^»(L.T^o) 

+ ((iVph) + l)G^<'»(w±wo) 



(A16) 



fact that one obtains 8 equivalent diagrams of the 7q 
order. 

Taking the steady state limit and after Fourier trans- 
formation, the different Keldysh components of the 
double-exchange (DX) self-energy are given by 



e-vib 



7q f du dv 



C3,C4 (A19) 



3. Normalisation of the vibron propagators 

The expression for the self-energy 5]£^^(ri,T2) can 
actually be recast in a Fock-like diagram with a renor- 
malised vibron propagator Dqt^^: 

Sf.vib(n, ra) = i7o' 2?o(ti, T2) G(ri, T2) , (A20) 
where the dressed vibron Vq{t,t') is given by 

2?o(t, t') = 72 / dT2dT3 Do{t, ra) n(r3, T2) i?o(T2, r') 

JCk 

(A21) 

and the polarization n(r, r') is given by the electron-hole 
bubble diagram: 



n(r,r') = -i G(r,T') G(r',r) 



(A22) 



2. Second-order self-energies 

The two- vibron process self-energies correspond to the 
two diagrams shown in Figure |3l i.e. the so-called 
double exchange (DX) diagram and GW-like (dressed 
phonon/ vibron DPH) diagram in which the vibron prop- 
agator is renormalized by a single electron-hole bubble 
polarization. 

The expressions for the self-energies for these diagrams 
are: 

4 r 

sf^b(n,T2) = / dT3dnG{Ti,T3)Do{n,T4)x 

JCk 

G{t3,T4)Do{t3, T2)G{t4, T2), 

(A17) 



and 



^0. 



S^vib(n,r2)^+^G(ri,r2)x 



/ dT3dT4Dn (ti , T3)G(t3 , T4)G(t4, T3)Do (T3, T2). 

JCk 



So, in principle, we already have all the ingredients to 
perform calculations using the fully dressed vibron prop- 
agator 2?(a;)^^-'^^ whose advanced and retarded compo- 
nents are given by 



(A23) 



(A18) 



The vibron spectrum is renormalized by the polariza- 
tion. Ther is a shift of the vibron energy related to 
^Q^ell^^"'^ and a finite linewidth of the peaks related 
to 7g3mn'f'"^°'^, instead of Dirac S peaks as in Eq. (|Bip . 

We have not performed the full vibron renormaliza- 
tion in the present work, because we want to study the 
effects of different diagrams of the same order, and order 
by order. I?(w) correspond to a partial resummation of 
all the bubble diagrams and contains all (even) orders of 
the electron-vibron coupling. The effect of this GW-like 
electron-vibron interaction will be presented in a forth- 
coming paper. 



Appendix B: Details of numerical calculations 



The factor 1/3 comes from the series expansion of the 
exponential in the time evolution operator Sck ^^id the 



In this section we provide more details about how cal- 
culations for the self-energies are performed. 
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FIG. 14: Equilibrium spectral functions for the off-resonant 
regime calculated with MCIST and NEGF-SCBA. MCIST 
calculations are performed at the Hartree-Fock level (BA). 
NEGF are done using the limit rj ^ 0^ for the imaginary 
part of the bare vibron GF and also by using a tiny, but 
finite value for rj (0.001 and 0.025). See main text for details. 
The other parameters are eo ~ +0.5 (electron transport), 
7o — 0.20 (upper panel), ujo = 0.40, toL = ton — 0.15. For 
a small enough value of r], one recovers the exact MCIST 
as well as the 77 ^ 0"*" NEGF results obtained from Hilbert 
transform (HT). Larger values of widen the peaks; however 
the spectroscopic information remains. For such values of r;, 
one can work with fewer grid points and makes second-order 
self-energy calculations tractable. 



First, we start by using the following decomposition of 
the retarded and advanced vibron propagator Dq"^: 



1 



1 



= P.V. — ^- P.V. ^ 



(Bl) 



and the definition of the Hilbert transformation (H.T.): 



■H{y) = H.T.[/(x)] = -P.V. r dx 

J-00 



x-y' 



(B2) 



One can then rewrite the Fock retarded (advanced) 
self-energy in terms of retarded (advanced) and lesser 
electron Green's functions only: 



+ o i'^lG< + G^]{cj + ^0) - H[G< + G'-]{lo - loo)) 



(B3) 



and 



SF.a / \ A 
e-vib(^) =70 



(iVph)G"(c^ - c^o) + ((A^ph) + l)G"(c.' + ^0) - 2 (G<(^ - ^0) - G%i^ - c^o) - G<(c<. + c^o) + G''(c<. + c^o)) 



('H[G< - GT(w + wo) - H[G< - G-^Kcj - wo)) 



r 



(B4) 



These expressions are in a form convenient for com- 
putation, since, as for S^ '^j,^ , they involve the direct 
evaluation of G"''°'^ on a energy grid shifted by ±wo in 
addition with the Hilbert transform of and C'" on a 
shifted energy grid. 

The Hilbert transformation is actually a conventional 
convolution product of the trial function with an appro- 
priate (l/a;)-like kernel. Calculation of convolution prod- 
ucts can be made faster, instead of scaling as the square 
of the number of grid points, by the use of FFT routines. 



However since FFT routines introduce artifical periodic 
boundary conditions, and one has to make sure that the 
functions, which decay slowly as l/w (the kernel of the 
Hilbert transform, and real part of G^'"'), have sufhciently 
small values at the boundaries of the energy grid. Oth- 
erwise, the corresponding discontinuities introduce spu- 
rious oscillations in the FFT. This means that one has 
to work with a wide energy grid and with a correspond- 
ingly large number of grid points to keep a good enough 
resolution. 
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Working with a large number of grid points iVpts > 10^ 
is not really a problem when calculating convolution 
products with FFT, since the number of operations re- 
duces from N^^^ to Npts log2 -/Vpts, or even when without 
using FFT routines. However, it is not possible to re- 
duce the calculation of the second-order DX self-energy 
to some sort of convolution product, and the number 
of operations then scales as N^^^. Such a scaling starts 
to make calculations seriously impractical when working 
with a large number of grid points. 

Hence, we have adopted another stategy which consists 
of introducing, in the definition of i'g'", a tiny but finite 
imaginary part iry instead of taking the limit i] 0+. 
This is just a natural step towards working with a fully 
dressed (renormalized) vibron propagator as already dis- 
cussed in Section IX3l 

By using this tiny but finite imaginary part iry, we avoid 
having to deal numerically with the Dirac (5-function and 
consequently with the Hilbert transform and all the prob- 
lems associated with the slowly decaying kernel. We can 
then work with a less wide energy grid and fewer points 
while keeping a good energy resolution. Hence the calcu- 
lations of the second order DX self-energy become more 
tractable. Of course, there is a price to pay for that: 
there is a lowest bound for the possible values of 77, and 
this lowest bound is strongly linked to the value of the 
grid spacing. By trial and error, we have found that rj 
must be at least equal to 2 or 3 times the grid spacing. 

Now we have to check and compare the results for 
the spectral functions obtained from these two differ- 
ent methods of calculation. Figure [14] show the spectral 
functions for the off-resonant regime at equilibrium and 
for weak/intermediate electron- vibron coupling. Calcula- 
tions have been performed with MCIST at the Hartree- 
Fock level (MCIST-SCBA) which serves as a reference 
calculation. The NEGF-SCBA calculations performed 
within the limit 77 — > 0^ by using Hilbert transform and 
FFT for the Fock self-energy (in the present case with a 
131072 grid points ranging from -30 to -1-30) are identical 
to the MCIST-SCBA results as expected. The NEGF- 
SCBA calculations performed with a tiny but finite imag- 
inary part irj have been done for different sets of param- 
eters. Only two are shown in Figure FHl 77 = 0.001 corre- 
sponding to a energy grid ranging from -10 to -1-10 and 
with 32769 points, and 77 — 0.025 corresponding to a en- 
ergy grid ranging from -10 to +10 and with 2049 points. 
The spectral function obtained with 77 = 0.001 is virtu- 
ally identical to the corresponding "exact" calculations, 
while the spectral function obtained with rj = 0.025 has 
slightly broadened peaks with a reduced amplitude, even 
though such a calculation has been performed with only 
2049 grid points — which is in the range of good values of 
grid points to make calculation of the second-order DX 
self-energy numerically tractable. In any case, the most 
important point is that the spectral information (peak 
positions) are not dependent on the finite value of rj. 

Now we also have to verify the infiuence of 77 on another 
transport regime for which we do not have any exact 
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FIG. 15: (Color online). Equilibrium spectral functions for 
the resonant regime obtained within NEGF (Fock self-energy 
only). Calculations are performed for different values of the 
imaginary part 77 (0.03, 0.01875, 0.00938, 0.001875, 0.000938) 
and for different number of energy grid points (2049, 32769, 
65537, 131073). The energy grid is ranging from -10 to -f 10. 
The inset shows a zoom around = of the rescaled spectral 
functions, giving then the linear conductance Gun at = 
0. The other parameters are eo = 0.0, 70 = 0.23 , cjq = 
0.40, ioL = %R ~ 0.15. Peaks in the spectral functions narrow 
with decreasing values of as expected. And one recovers a 
perfect linear conductance Gun — Go s,t T — 0. 



reference calculations, but for which we may expect to 
recover conventional properties of an interacting Fermi 
liquid at very low temperatures. 

Figure [TS] shows the equilibrium spectral functions for 
the resonant regime. Calculations are performed with 
only the Fock self-energy and for different values of 77 
(0.03, 0.01875, 0.00938, 0.001875, 0.000938) and for dif- 
ferent number of energy grid points (2049, 32769, 65537, 
131073). The energy grid is ranges from -10 to -1-10. As 
expected, the main central peak and the vibron side-band 
peak narrow with increasing amplitude when the value of 
77 decreases, until reaching an asymptotic behaviour for 
77 ~ 10~3. 

The most interesting result is given by the behaviour of 
spectral function at the equilibrium Fermi level fi"'^ = 0. 
Such a behaviour is exemplified in the inset of Figure 
[T5] The inset shows the spectral function rescaled by 
2r{fi'"^) so that the value at w = corresponds to the 
linear conductance Gun (see Section fill Ep . One can see 
that Giin goes towards the unit of conductance as 77 — > 0. 
Such a limit corresponds to a perfect conductance, as 
obtained for the corresponding non-interacting system. 

This result shows that even in the presence of electron- 
vibron interactions, the linear conductance is perfect, 
as one also obtains for the non-interacting case. The 
electron-vibron interaction (as least at the Fock level) 
does not renormalize the linear conductance. For the 
fully electron-hole symmetric system (no Hartree self- 
energy), we have checked that the value of QmT,l._^-y^{u}) 
is zero at a; = and at T = (at equilibrium). The 
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real part of the interacting self-energy 5ReEg_^jj^ = at 
a; = by definition because it is an odd function of 
Lo. Our electron- vibron interacting system conserves its 
Fermi-liquid properties at zero temperature and equilib- 
rium, as also obtained for interacting electron systems 
(see for example Refs.jssl-loU). We expect that such a 
behaviour also hold for other electron-vibron diagrams 
than the Fock diagram. Out of equilibrium, the system 
may lose its conventional Fermi-liquid properties, as we 
have already shown in the main part of the paper. 

In conclusion, we want to stress that by introducing 
a finite, but tiny, value for the imaginary part 77, one 
might not get the exact linear conductance value. How- 
ever this is a numerical artifact which can be controlled 
by increasing the number of energy grip points and re- 
ducing the value of 77. Nonetheless we expect that the 
relative effects of different electron-vibron diagram are 
correctly obtained from our calculations, whatever the 



small values of 77 are. And of course, in order to perform 
a correct analysis of such effects, one should only compare 
results for either spectral functions or linear conductance 
obtained numerically with the same value of 77. 
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